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This work studies the issue of quadrotor trajectory tracking control in presence 
of disturbances and model uncertainties. The paper starts by extracting the 
kinematics and dynamics models of the quadrotor. This results in the motion 
equations, which eventually serve as a blueprint for creating the suggested 
smart flight control scheme. Secondly, an enhanced backstepping controller 
(BSC) is developed and tested to keep the quadrotor tracking the desired 
trajectory both in steady state and in presence of disturbances. Finally, BSC 
beside two other controllers: sliding mode controller (SMC) and proportional 
derivative controller (PDC) are implemented in MATLAB/Simulink and the 
obtained results are compared and conclusions are extracted. Therefore, it is 
established that PDC is not robust to disturbances as noise will be amplified 
due to the derivative term. Whereas, although SMC is robust to parameter 
variations and disturbances; however, it is not continuous which may affect 
the actuators due to the increased gains which may saturate them. In contrast, 


BSC requires too many tuning parameters; however, it ensures Lyapunov 
Stability and does not depend on the system as it does not involve cancelling 
system nonlinearity. Moreover, BSC results are 10!” better than the results of 
the two other controllers. 
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1. INTRODUCTION 

Unmanned aerial vehicles drones (known as UAVs) are self-flying robots. In several circumstances, 
it is preferable that the system achieves significant movements, moves both vertically and horizontally, and be 
as tiny as feasible. With the advancement of UAVs’ technology, their applications are extended to surveillance, 
healthcare, agriculture, civil and military fields, transmission line inspection and energy delivery [1]-[3]. 
Drones provide a number of fundamental advantages over conventional aircraft systems, including increased 
flexibility, reduced price, fewer radar scans, a longer lifespan, and no harm to the operator's life. The quadrotor 
rolls and pitches in the direction of the slowly rotating motor. By dividing the thrust into two directions, 
represented by the roll and pitch angles, linear motion is generated. However, actually, the drone is considered 
a nonlinear system with resilient coupling and under-actuated characteristics [4], [5]. It can also be significantly 
impacted by outside disturbances like wind gusts. Therefore, the quadrotor provides a hard control problem 
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because of its very unsteady characteristics [6]. Hence, the need for an effective control method arises from 
the peculiar nature of the aerial vehicle. 

Recently, proportional integral derivative (PID) controller based on optimization algorithms [7]-[9] 
and fuzzy logic [10]-[13] have proved their efficiency in terms of controlling such nonlinear, unstable systems; 
however, their implementation on quadrotors was very limited because of the presence of external disturbances 
which are hard to be predicted. In the literature, several techniques were developed for vehicles’ control in 
general [14]-[16] and for quadrotors’ control in particular such as sliding mode controller (SMC) [17], SMC 
approach built on backstepping controller (BSC) which produces high performances and faster response [18], 
[19], sliding controller combined with state observer compensate for uncertain nonlinear components [20]. As 
well, Backstepping based nonlinear control and adaptive SMC [21], PID control [22]- [25], nonlinear adaptive 
and predictive controllers [26]-[28], plus the controller eliminating active disturbances (ADRC) which deals 
with the study of path following performances where considering outer factors [4]. Besides, adaptive 
proportional integral derivative (APID) controllers ensured faster convergence speed and minimum overshoot 
in the existence of constraint variations and outside instabilities [29], [30]. 

The BSC has proved its efficacy in the control of nonlinear systems compared to the other controllers, 
especially with the influence of exterior forces [5], [31]. Research by Fahmizal et al. [32] the three approaches 
(PID, SMC, and BSC) were applied to a standalone photovoltaic system's single-phase voltage inverter; the 
authors proved that BSC has the potential to perform better than the other two controllers. Hence, our work 
explores the issue of quadrotor trajectory tracking control with the influence of disturbances and model 
uncertainties. First, the kinematics and dynamics of the quadrotor are described in detail to aid in understanding 
the behavior of the drone. The motion equations that result from this [33] are then employed as a blueprint for 
creating the suggested smart flight control strategy. In this work, an easy-to-implement enhanced BSC is 
developed; then, tested beside two other controllers: sliding mode (BSC) and proportional derivative controller 
(PDC) to keep the quadrotor tracking the desired trajectory both in steady state and in the existence of outside 
instabilities. Finally, the three controllers' results are compared to determine which is the best. 


2. METHOD 

The distance between the body and the earth frames r = [x y z]? reveals the precise location of the 
quadrotor's center of mass. The roll (¢), pitch (0), and yaw (y) angles, which denote rotations about the X, Y, 
and Z axes, respectively, characterize the orientation of the quadrotor. The orientation of the quadrotor is 
expressed by the rotation R from the inertial to the body frame. Supposing that the rotation occurs in the 
following order: pitch, roll, and yaw; thus, R is assumed by (1) [34], [35]: 


cos(@) cos (y) cos(y)sin(8) sin($) — cos($) sin (y) cos($) sin(@) cosQp) + sin($) sin (y) 
R = |cos(0) sin (Y) sin($) sin(@) sin) + cos(@) cos (Y)  cos($) sin(8) sinQp) — sin(@) cos (Y) | (1) 
—sin (8) sin(@) cos (8) cos(@) cos (8) 
in order to connect the inertial frame-measured Euler rates h = [@ 6 w]” and the rates of the angular 
body w = [P 4 1]', the following conversion is required [34]:  — R,1j. With: 


R, = 0 1 sin(Q) (2) 


cos(0) 0 —cos(@) sin(@) 
bon 0  cos(Q)cos(0) 


small angle assumptions are made around the hover location where cos($) ~ cos(@) ~ 1 and sin($) ~ 
sin(0) « 0. R, can therefore be reduced to an identity matrix J [36]. 


2.1. Quadrotor dynamics 
Quadrotor's movement may be split into 2 control systems: a rotating sub-system (roll, pitch, and 


yaw) and a translating sub-system (x, y, and z coordinates). The revolving sub-system is completely actuated 
whilst translating sub-system is under-actuated [37]. Newton Euler equation can be formulated by (3): 


Lx 0 0 
Epe eel — seu? o 


With: F: net force acting on the quadrotor 
T: net torque 
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Ig3: 3 x 3 identity matrix 

14: moment of inertia 

m: the mass of the quadrotor 

a: the center of mass's linear acceleration 
a: angular acceleration 


2.1.1. Rotational subsystem 
In the body frame, based on the Newton-Euler method the rotational equations of motion are derived 
using the general formalism (4): 


Iù +w x Iw = Mg (4) 


Where: T: is the inertia matrix of the quadrotor 
w: is the angular body rate 
Msg: are all the moments acting on the quadrotor in the body frame 
Io and w X Io: represent the rate of change of angular momentum in the body frame 


2.1.2. Matrix of inertia 

The quadrotor's inertia matrix is a diagonal matrix, and because of the symmetry of the quadrotor, the 
off-diagonal elements, which make up the inertia product, are zero as shown in (3). Ixx, Iyy, and Izz are the 
body frame's primary axes' respective area moments of inertia. Before defining the last term (Mp) of (4), two 
physical effects have to be defined: the rotor-generated moments and aerodynamic forces. Each rotor produces 
a force called the aerodynamic force F; and a moment called the aerodynamic moment Mi as a result of rotation 
[32]. They are given by (5) and (6): 


Fy = 5 pACpr?02 = K0} (5) 


M; = 5 pACpr?.0? = Ky? (6) 
Where: p:is the air density 

A: blade area 

Cp, Cr: aerodynamic coefficients 

r: radius of blade 

Q;: angular velocity of rotor i 

Ky, Ky: are the aerodynamic force and moment constants respectively. 
Figure 1 depicts the forces and moments acting on each of the quadrotor's four rotors. 


Figure 1. Moments and forces affecting the quadrotor’s motion 
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It is obvious that each rotor creates a moment M; with a direction opposite to the corresponding rotor's 
directions and an upward thrust force Fj. By relating the moments around the axes of the body, the right-hand 
rule combined with the body frame's axes facilitate the deriving of the equations of Mg [35], which are written: 


IK, (—23 + 02) 
Ms = IK, C02 — 02) (7) 
Ky(—0? + 02 — 03 03) 


so, by replacing Mz in (4) and from (1) we can get: 


$ cos(@) 0 sin(@) p 

,| _ |sin(@) tan(@) 1 -cos(0)tan(Q) 

0 = __sin(@) 0 cos(@) : A 
y cos() cos($) 


2.1.3. Motion's translational equations 
Using the Earth's frame and Newton's second rule of motion, the linear equations of motion are given 


0 
»-| 0 


—mg 


by: 


=E + RF; (9) 


Where: r=[X Y  Z] define the distance of the quadrotor from the inertial frame 
g: acceleration due to gravity g = 9.81 m/s? 
F,: drag forces 
Fp: nongravitational forces. To translate the thrust forces from the body frame into the inertial frame, 
non-gravitational forces operating on the quadrotor are multiplied by the rotation matrix 
The nongravitational forces Fp, the resisting force F, and the drag moment M, could be stated as: 


0 
Fg = 0 (10) 
Ky(0? 02 0$ + 01) 
F, = Kt (11) 
M, = K,n (12) 


Where: K;: a constant matrix called the aerodynamic translation coefficient matrix 

T: the derivative of the position vector r (velocity of the quadrotor) 

K,: a constant matrix called the aerodynamic rotation coefficient matrix 

7: the Euler rate 
From (4) can be written as: 

Io + oxlo-ctM;—-Mg-M, (13) 
2.2. State space representation 

Due to a lack of sensors, state variable measuring is expensive and challenging in reality [38], [39]. 
Thus, quadrotor's state vector is defined as: 


X = [x1 X2 Xs X4 Xs X6 X7 Xg Xo X10 X11 X12]. (14) 


this is translated into the degrees of freedom of the quadrotor presented in (2). The quadrotor's position in 
space, together with its angular and linear velocities, are all indicated by the state vector as: 


. . . 2 A 1T 
X=|ġpġ0ġypýzżxžyy] (15) 
having: U, = K; (Q2 + 03 + 02 + 02), U, = K-(—03 + 03), U; = K; (02 — 02), 
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U, = Ky(-02 + 02 — 02 + 02) 


-I L,-I L.—I l l l 
yy ZZ ZZ XX XX yy 

=] = ————————, = —— _,b = —,b = —,b = — 
i : I a lyy ye I ; Trex - lyy i Iz 
X, = $ = X2 
X, = Q = a1 X4X6 — a5x4Q + DU, 
X; = 6 = x4 
X4 = 0 = X2Xg03 + 4X20. + b2U3 
Xs = Y = X6 
Xg = Y = X2x40s + b3U, 
ky =Ż= Xg (16) 


P T 1 
Xg = Z— —g + —cos(x4) cos(x;) 
m 
Xo = x = X10 
j^ Us ; 
Xio = X = — (sin x, sin xg + cos x, sin x4 cos xg) 
m 
X11 — Y = X312 


. T 1 . . . 
a ae (— sin x4 cos xg + cos x, sin x3 sin Xs) 


2.3. Quadrotor's backstepping control 

Backstepping control is a recursive Lyapunov based control technique for systems in strict feedback 
form, it works by cascading a number of sub-systems into the main system. The control rules are then 
developed, one for the entire system and one for each subsystem, in decreasing order. The elaborated rules are 
given in the following sections. 


2.3.1. Roll control 
Taking into account the first stated sub-system: 
ži =x% 
m = a4X4Xg — a5x4CQ + bU; CN 
Step 1 
The following is an expression for the error £, between the desired and actual roll angles: 
& =x] X (18) 


Consider the Lyapunov function 


1 
Vi = zë 


Therefore, V is derived through x; variable, V, is calculated as follow: 
V, = & E 
With: 
E Ni cde ea 
Choosing 
& —-K48 
(where: K4£4 positive definite function), we get: 
ae. TES 


Step 2 


Denoting €, the error among the actual roll angle rate and the desired one, thus: 
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Defining: 


V = V; + d 
2 
being a potential Lyapunov function 
V, = Vi + £28, = £161 PES = —K] E2 + £, (£1 + žÎ — x4X601 + x4Qa, — b, U) 
Supposing 
E (£1 + X2 — x4xga4 + x40a; — bU) = -K3€2 


1 
U, = p, (Kee + Ey — X4X$04, + a5x40 — kıx3) (19) 


2.3.2. Pitch control 0 
Taking into account the second stated sub-system: 


{ fa He (20) 


X4 = X2X6a3 + 4X20. + b2U3 


Step 1 


Considering €, the error among actual pitch angle rate and desired one and which can be obtained via: 
£z = XÍ — X3 > é c e 
Utilizing Lyapunov stability, pick: 
V (e3) = La 
2 
The system trajectory is ensured to check the following condition if V is negative: 


Y, = £363 = ez (x2 - X4) « 0 


Then: 
qox = —Kz£3 > > xf = ži auN e. 
Step 2 
The error: 
PE E r P E E E a T 
V, = V; + 28 > V, = &383 + £484 = & (X2 — x4 + £4) + &4(£4) 
Putting: 
£4(£4 — X2Xga3 — 44X2 + b2U3) = —k4e2 
U3 = s; Chats — £3 + X2X6ü3 + a4x2 0 + kax4) (21) 


2.3.3. Yaw angle control v 
Now think about the third subsystem listed beneath: 
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Xs = X6 
P = X2X4đ5 + b4U, (22) 


Step 1 
£g describes the error among the actual yaw angle rate and the desired one. Accordingly: 


£s = xf — X; > £; He — ý; — 


with, Lyapunov function is 


V(és) = Led >> V; mel 


thus; 
aoc = —K;£; > x = ee + E 
Step 2 
The error 
E6 = XÊ — X6 > Xe = XË — £6 
Ve = V; + Lez => Vs = Ezé + Ecés > Ve = Es (tt — xf + E6) + S 
Naming: 
E6 (E5 — X2X40g — kgxg — b4U4) = —k6£6 
U, = 5 (kets + £s — X2X45 — ksXg) (23) 


2.3.4. Altitude control 
In (24) represents the fourth subsystem: 


X; = Xg 
ls = COS(X1) cos(x2) U, ~g (24) 


m 


Step 1 
€- describes the diference amongst the actual position z and desired one: 


Eo Sm ey SSO UE e 


Lyapunov function is 


V(e7) = Lej >> on 


Then, 

Xo x = —K;£&; > > xf um WE E 
Step 2 

The error 

Eg = XË — Xg Ə Xg = Xİ — £g Ə Éz XS — kg. 
Thus, 
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1 n 
Vs = V, + P => Vg = £;&; + EgEs 
r d d H d U, 
Ez; = €7(x4 — xf) + €7&g + Eg(Eg) = £g | €; + X8 | g — 08 63)cos(x) 


Putting: 
U. 
£g (s - xd - (o - = cos(a) co 62) = — kgeg 


Thus: 


U, = —À—— (-€, + g — kg £g kx; ) 


cos(x4)cos (x2) 


2.3.5. Control y position 
Formula (26) denotes the fifth sub-system: 


X9 = X10 
: Ui 
X10 — Uy a 


Step 1 


Defining £, the difference concerning the actual position y and desired one: 


Pe — Xo > é = Ži — žo 


The Lyapunov function is 


1 ; . 
V(&) = 56 > Vs = €9€&5 


Then, 
xê — Xio ee 6 > > xÅ = xf + NES 
Step 2 
Then, the error is: 
mue IA — 2d 4 
E10 = Xio — X10 > £10 = X10 T X10 
V S V 1 2 V = H H H —-— dá d . 
10 = Vo + 2510 > V = £gég + E10€19 > > — &gég = £g(X$ — X19 + £10) + &10(£10) 
U 
„d 1 _ 2 
£10 e + X79 + — Uy) = -K79€{o 
m 
Therefore: 


m 
Uy = Ui (-€9 — k10 &1o- K9X40) 


2.3.6. Control x position 
Underneath, (28) denotes the final sub-system: 


Xii = X12 

TELA 

X12 — Un 
Step 1 
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(26) 


Q7) 


(28) 
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Defining € , the error between the actual position x and the desired one: 
&1 = xf -xu > ài -ih-iu 

Lyapunov function is 


1 , l 
V(&) = 5 et > Vg = £g 


Thus, 
„d = d _ d 
X4 — X12 = —Kyy 611 > xi? = X5 + Ky 11 
Step 2 
Set the error £4; as: 
d. E sg —.d x 
€12 — X12 — X12 €12 — X12 — X12 
Thus, 
NM 
V = Vi + 7fi2 > Vi? = £11811 + £12612 
Then 
U. 
„d 1 = 2 
€12 (ss + xiz + -Lu,) = -Ky2€12 
m 
Therefore: 


U, = D (7&1 — kaiz E12 7 K11X12 ) (29) 


3. RESULTS AND DISCUSSION 

MATLAB/Simulink is the most useful software used to test the behaviour of nonlinear systems and 
to validate the results of the recently developed controllers [38], [40], [41]. Therefore, model simulation is 
done using MATLAB/Simulink program to verify our suggested control approach. To achieve that, a quadrotor 
model is designed and controlled by three different controllers which are: PD, sliding mode and the BSC. 
Hence, each one of them is tested with and without disturbances to track a quadrotor's trajectory of radians 
R=8 shown in Figure 2. The chosen quadcopter's mass m=200g; whereas, the area moments of inertia /xx, Iyy 
and Izz are 0.00025, 0.000232 and 0.0003738 respectively. The quadrotor has to track the trajectory defined 
by the time functions: X = 8 sin (0.1t), Y = 8 cos (0.1t ) and Z = 0.2t. 


3.1. Without disturbance 

Figure 2(a) is an illustration of the trajectory tracked by the quadrotor. The next figures show the 
position, orientation (see Figures 2(b)-(d)), trajectory errors (illustrated in Figures 3(a)-(c)) and control inputs 
in the absence of disturbances generated by the three controllers: PDC (illustrated in Figure 4(a)), SMC (as 
shown in Figure 4(b)) and BSC (see Figure 4(c)) (see in appendix). Discussion: without any disturbance, the 
steady state error from the SMC was the smallest followed by PDC then the BSC for the displacement across 
both the X and Y axes. In the event that no disturbances, SMC has proven its efficiency even in the previous 
work [4], [20], [42]. However, for the displacement on the Z axis, the SMC did 10? times worse than the BSC, 
with the PDC doing the best here. And finally, for the steady state error for yaw angle, the BSC also was 10!” 
better than the PDC, and the SMC and being lasting this category. 


3.2. With disturbance 

The next figures show the position and orientation (see Figures 5(a)-(c)) (see in appendix), the 
trajectory error Figures 6(a)-(c) (see in appendix) and control inputs (shown in Figures 7(a)-(c)) in the presence 
of a ramp disturbance of the vector F=9t i+9t j+ 9t k (N), starting from the 10" second of the simulation with 
a force limit of 9 N. Simulation are generated by applying the two controllers: PDC (shown Figure 5(a) and 
Figure 7(a)) and SMC (shown in Figure 5(b) and Figure 7(b)) versus BSC (illustrated in Figure 5(c) and 
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Figure 7(c)) (see in appendix). Discussion: now, in the case of disturbance. For small disturbance, the steady 
state error for every controller kept the same order as in case of no disturbance. But as the disturbance increased 
the controllers could not keep the quadrotor in trajectory anymore. The first controller that collapsed was the 
SMC. Making the PDC the best controller of the three under disturbance, but even this last one also collapsed 
after adding 133.33% of the first disturbance force. Then, BSC proves that it is the best to resist for all kind of 
disturbances; but, after adding about 50% of the previous disturbance force, the BSC collapsed too. However, 
the simulation result obtained from the implementation of the proposed BSC is very satisfactory compared to 


the previous works which show some complexity of the analytical inference and a considerable dynamic error, 
especially in presence of noise [4], [5], [43]. 


x(m) va time(s) roll (in deg) vs time(s) 


0 2 4 & 8 10 


Uo 23 4 © 9 10 
zim) vs time(s) 


>» © © 9$ 
xt [in dog) vs time(s) 


Figure 2. Attitude, heading and position reference (a) trajectory tracked by the quadrotor, (b) PDC, (c) SMC, 
and (d) BSC 
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0.02 0 -0.02 
0 20 4 60 8 10 0 2 4 6 8 100 0 20 4 6 8 10 0 20 4 ü 8 100 
» error in z(m) vs time(s) error in psi(yaw)(deg) vs time(s) error in z(m) vs time(s) error in psi(yaw)(deg) vs time(s) 
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-0.02 
5 -0.01 
-0.04 
ot 0.02 -0.06 
0 20 4 60 8 10 0 20 4 60 8 100 0 20 40 60 H8 
(a) 
error in x(m) vs time(s) T Roll error vs time(s) 
02 1 | 
0 J 0 
pp ——— ————— "n —— ——— —— 
0 20 40 60 80 100 0 20 40 60 80 100 
error in y(m) vs time(s) Pitch error vs time(s) 
02 
02 
| YN /N 
-0.2 + + -0.2 n n 
0 20 40 60 80 100 0 20 40 60 80 100 
ñ error in z(m) vs time(s) emt psi(yaw)(deg) vs time(s) 
0 
5 2 
E 
0 
0 20 40 60 80 100 0 20 40 60 80 100 
(c) 


Figure 3. Altitude, heading, and position error signals (a) PDC, (b) SMC, and (c) BSC 


4. CONCLUSION 

This paper proposes a controller for quadrotor's trajectory tracking. First, kinematics and dynamics 
modeling of the quadrotor are developed which yield the motion's equations. Basing on this, an enhanced 
backstepping scheme is developed for position and orientation subsystems where stability analysis is ensured 
by Lyapunov concept. 

It is established that the position, orientation and attitude path following errors can rapidly converge 
to slight values with all controllers. In case of non-external disturbance, BSC shows good control of the yaw 
angle and the altitude of the quadrotor comparing to the two other controllers (SMC and PDC). Moreover, in 
the case of presence of disturbances, for small disturbance, each controller's steady state error maintained the 
same order as in the absence of any disturbance. However, as the disturbance increase the controllers could not 
keep the quadrotor in trajectory anymore. Numerical and simulation results confirm that BSC is the last one 
that collapsed which confirm the robustness and efficacy of our constructed enhanced control strategy. Our 
perspectives involve improving the proposed approach; then, test and design it on real experimental quadrotor 
platform. 
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APPENDIX 
à U1 vs Time(s) «10%  U2vs Time(s) 25 U1 vs Time(s) 4X 10?  U2vs Time(s) 
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Figure 4. U1, U2, U3 and U4 vs. time, (a) PDC, (b) SMC, and (c) BSC 


x(m) vs time(s) roll (in deg) vs time(s) 


Uo 2 4 & 8 100 0 20 4) & © "0 0 23 4 © 95 10 
yim) vs time(s) pitch (in deg) vs time(s) yim} vs time(s) 


UD 2 a & 8 w — 0 2 æ & æ w X2 4 & B 10 
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Figure 5. Altitude, heading and position reference measurement vs. actual measurement (a) PDC, (b) SMC 
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Figure 6. Altitude, heading and position error signals (a) PDC, (b) SMC, and (c) BSC 
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Figure 7. U1, U2, U3, and U4 vs. time (a) PDC, (b) SMC, and (c) BSC 
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